Non-Gaussian dynamics from a simulation of a short peptide: 
Loop closure rates and effective diffusion coefficients 
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Intrachain contact formation rates, fundamental to the dynamics of biopolymer self-organization 
such as protein folding, can be monitored in the laboratory through fluorescence quenching mea- 
surements. The common approximations for the intrachain contact rate given by the theory of 
Szabo, Schulten, and Schulten (SSS) [J. Chem. Phys. 72, 4350 (1980)] and Wilemski-Fixman (WF) 
[J. Chem. Phys. 60,878 (1973)] are shown to be complementary variational bounds: the SSS and 
WF approximations are lower and upper bounds, respectively, on the mean first contact times. As 
reported in the literature, the SSS approximation requires an effective diffusion coefficient 10 to 100 
times smaller than expected to fit experimentally measured quenching rates. An all atom molecular 
dynamics simulation of an eleven residue peptide sequence in explicit water is analyzed to investi- 
gate the source of this surprising parameter value. The simulated diffusion limited contact time is 
approximately 6 ns for a reaction radius of 4 A for solvent viscosity corresponding to that of water 
at 298 K and 1 atm (rj = 1.0 cP). In analytical work, the polymer is typically modeled by a Gaussian 
chain of effective monomers. Compared to Gaussian dynamics, the simulated end-to-end distance 
autocorrelation has a much slower relaxation. The long time behavior of the distance autocorrela- 
tion function can be approximated by a Gaussian model in which the monomer diffusion coefficient 
Do is reduced to Do/6. This value of the diffusion coefficient brings the mean end-to-end contact 
time from analytical approximations and simulation into agreement in the sense that the SSS and 
WF approximations bracket the simulated mean first contact time. Attention to the non-Gaussian 
nature of the dynamics has direct implications for the development of improved analytical models. 
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Introduction 



Intrachain contact formation is fundamental to our un- 
derstanding of biopolymer dynamics involved in molec- 
ular self-organization. Motivated primarily by the rel- 
evance to protein folding, several groups have recently 
measured the quenching rates of the terminal residues on 
short peptides, using fluorescent probes with short range 
quenching. una In these experiments, the fluorescence de- 
cay rate characterizes the rate of contact formation in 
peptides for which the dynamics should be less compli- 
cated (and hence, easier to interpret) than the folding 
dynamics of natural protein sequences. Lapidus et ala 
use a simple polymer model with a single dynamical vari- 
able to understand the measured rates, finding that the 
effective diffusion constant needed to fit the data is ap- 
proximately 30 times smaller than the expected relative 
coefficient for free diffusion of the two terminal residues 
(assuming this to be 1.5 x 10 -5 cm 2 /s). A small effec- 
tive diffusion coefficient appears to be a rather general 
result, consistent over different molecules and types of 
probes. Many years ago, Haas et ala also reported sur- 
prisingly small diffusion coefficients using similar poly- 
mer models (see also Ref. ^|), but here the quenching 
mechanism of the terminal monomers of the peptides is 
through the longer range fluorescence resonance energy 
transfer (FRET). Recently, Wallace et. ala applied this 
analysis to FRET quenching experiments in DNA, con- 
cluding that the effective monomer diffusion coefficient 
was 1000 times slower than reported for peptides. 



Interpretation of this small effective diffusion coeffi- 
cient is particularly relevant to modern protein fold- 
ing research. The energy landscape theory of protein 
foldingO along with advances in experimental techniques 
that probe faster timescales has encouraged some re- 
searchers to focus on proteins with tho-.|Simplest ki- 
netics, the so-called fast folding proteinsflu There has 
been considerable progress in understanding these fast 
folding proteins from both experimental and theoreti- 
cal approaches. For this group of proteins, the transi- 
tion state ensemble probed experimentally through mu- 
tational studies pioneered by FershtJiS can be modeled 
fairly accurately by simple polymer based models with 
interactions between distant .residues specified by the na- 
tive state topologyJBll3't3oEa Success of modeling the 
barrier crossing dynamics that determine folding rates 
is far less clear; thy fp w th e ories of the folding dynam- 
ics proposed so fart§0£3o could greatly benefit from 
more experimentally supported microscopic parameters 
characterizing the peptide dynamics. Such fundamen- 
tal timescales and model parameters can be obtained in 
principle from measured intrachain quenching rates and 
more accurate theoretical modeling. 

Effective parameters are relevant to specific models. A 
common model for polymer dynamics approximates the 
polymer backbone connectivity by harmonic bonds be=. 
tween effective monomers. As emphasized by Zwanzig,E2l 
a general harmonic potential energy can be defined 
through the static correlations between the monomers. 
Interactions between distant monomers do not enter ex- 
plicitly into the model except possibly (and indirectly) 
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through the static correlations defining the Gaussian po- 
tential. The dynamics of these phantom chain models 
obey Gaussian statistics. In particular, the dynamics of 
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the end-to-end vector R(t) is completely described by the 
Gaussian Green's function (assuming isotropy) 
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Here, <p(t) is the equilibrium vector correlation 

(R(t) ■ R(0)) 



(R 2 



(2) 



and {R 2 } is the equilibrium mean square end-to-end dis- 
tance that characterizes the size of the polymer through 
the Gaussian equilibrium distribution 



P cq (R) 



2ir (R 2 



7 exp [-/?[/(#)], 0U(R) = ^- 

(3) 

As long as the Gaussian nature of the dynamics is 
retained, the only relevant difference between different 
Gaussian polymer models is through the static size (R 2 ) 
and the equilibrium pair correlation <fi{t). Neglecting fric- 
tion with memory,E3EII <f>{t) is typically expressed as a 
sum of exponentials with a wide distribution of decay 
times. Another common approximation treats R{t) as 
the only dynamical variable, replacing cj>(t) in Eq. (Q) with 
a single exponential. This single variable approximation 
to the multi-bead dynamics is the one that has been used 
to analyze the fluorescence quenching experiments. 

Some explanations for the small effective diffusion co- 
efficient required to fit measured quenching rates can be 
incorporated into Gaussian models. For example, Haas 
et am suggested that the reduced effective diffusion con- 
stant may be evidence for internal friction of the chain 
(independent of the solvent viscosity). Internal friction 
designed to capture the timescale of dihedral angle basin 
hopping has been introduced through a modified friction 
matrix while— retaining the Gaussian description of the 
dynamics .E3'E3 Alternatively, Lapidus et aln noted that 
the reduced effective diffusion constant is reminiscent of 
diffusion in a one-dime«sional rugged potential well as 
discussed by Zwanzig.C-3 Within the single variable ap- 
proximation, averaging over the rugged energy poten- 
tial can lead to dynamics described by diffusion on a 
smooth potential with a reduced effective diffusion con- 
stant. Both local chain dynamics or energetic ruggedness 
can provide additional friction that effectively slows the 
dynamics of the Gaussian chain. A main focus of this 
paper is to determine whether a Gaussian chain with 
modified parameters can capture timescale of the poly- 
mer dynamics suggested by the small effective diffusion 
coefficient. 



There are certainly other mechanisms that could slow 
the dynamics, some even leading to glassy behavior. The 
ruggedness in Zwanzig's one dimensional model are local 
extrema in the potential along the reaction coordinate, 
sometimes referred to as "longitudinal ruggedness." Al- 
ternatively, the Bryngelson and Wolynes theory of pro- 
tein folding dynamicsES accounts for "transverse rugged- 
ness" (the coupling to hidden or additional degrees of 
freedom) that may induce kinetic traps along the reac- 
tion coordinate. Interactions between monomers along 
the polymer chain can also generate a form of friction. 
Mode-coupling theories predict that the dynam^csxif ran- 



dom heteropolymers may become non-ergodic,l sim- 
ilar to spin glass dynapics below the dynamical glass 
transition temperatureO Although one would not ex- 
pect a short peptide to exhibit glassy behavior, similar 
mechanisms may still be relevant to some degree in this 
case as well. 

The general problem of intrachain quenching is a dif- 
ficult one that has not been solved explicitly even for a 
Gaussian chainH3 Two widely used approximations are 
the cloptre approximation proposed by Wilemski and 
FixmanEl (WF) and the single variable-approximation 
solved by Szabo, Schulten, and SchultenEl (SSS). In this 
paper, it is shown that the diffusion limited reaction rate 
from these two approximations are related through a vari- 
ational expression: the WF and SSS approximations are 
upper and lower bounds, respectively, to the mean first 
passage time of end-to-end contact of the original chain. 
Since the single variable approximation used to interpret 
the experimental measurements is a lower bound, it is 
possible that the SSS approximation requires a small ef- 
fective diffusion constant to fit the experimental quench- 
ing time because it is a poor approximation to the exact 
mean contact time of the multi-bead chain. Similarly, the 
WF approximation (had it been used) requires a larger 
diffusion constant to fit the data, since it gives a longer 
contact time than the SSS approximation. 

The accuracy of these approximations have been as- 
sessed by the comparison to simulations of a Rouse 
chain (raopomers sequentially connected by harmonic 
springs) The results reported by Pastor et a/.c3 

show that the SSS approximation underestimates the 
Rouse chain mean contact time by roughly a factor of 
three (although this factor is just an example as it de- 
pends strongly on chain length and contact distance). 
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Consequently, the accuracy of the lower bound may at 
least account for part of the small effective diffusion con- 
stant that needs clarification. 

In order to develop improved models informed by mea- 
sured quenching rates, it is crucial to clarify the way 
in which the effective diffusion coefficient reflects the 
timescale of the polymer dynamics. This question is 
addressed in this paper through an analysis molecular 
dynamics simulation of a realistic peptide. To connect 
with analytic polymer models, a general Gaussian chain 
is chosen to define the effective diffusion coefficient and 
compare against the simulated dynamics. 

The peptide sequence C(AGQ) 3 W (denoted by CW 3 ) 
studied experimentally by Lapidus et allu is simulated 
using all atom molecular dynamics in explicit solvent. A 
closely related work by Yeh and Hummero was recently 
published that presents extensive all atom simulations of 
the shorter peptides C(AGQ) 1 W and C(AGQ) 2 W. The 
simulated mean quenching time for CW3 is the same 
order of .magnitude as found previously for the shorter 
peptides:E3 for diffusion limited quenching with a reac- 
tion radius of 4 A we find mean contact time of approxi- 
mately 6 ns (adjusted to solvent viscosity ?y wat = 1.0 cP). 
The experimentally measured quenching time for this 
peptide (on the order of 100 ns) is many times longer 
than the diffusion limited contact time,EJ suggesting that 
the quenching rate for this system is reaction limited 
not diffuaipn limited. Indeed, the analysis of Yeh and 
HummerEj shows that the simulation is consistent with 
experiment assuming short distance quenching with a fi- 
nite rate.EZI While this is very important for interpreting 
these particular experiments, the fundamental time scale 
that one wishes to ultimately understand in the folding 
context is the polymer dynamics of intrachain contacts, 
not the quenching properties of the probes. The simu- 
lation of CW3 presented in this paper as well as those 
presented by Yeh and HummerE3 show that a reduced ef- 
fective diffusion coefficient is still required to fit even the 
shorter diffusion limited quenching time. This is the fo- 
cus of the present paper as the reaction limited quenching 
is a rather separate issue. 

The analysis of the present simulation suggests that 
the small effective diffusion coefficient does not origi- 
nate from local internal friction or longitudinal rugged- 
ness that can still be modeled in an otherwise Gaussian 
formalism. The accuracy of the single variable approxi- 
mation to the multi-bead quenching time does not seem 
essential to the discrepancy between calculated and mea- 
sured rates either. For this simple reaction, the prima 
facie reaction coordinate of the quenching would be ex- 
pected to be the end-to-end distance R(t) — \R(t)\. Con- 
sequently, R(t) is also the relevant quantity requiring an 
accurate dynamical description. Due to Eq.(|I]), the pair 
correlation (R(t)R(0)) for a Gaussian model can be found 
once (i? 2 ) and (f>(t) are specified. From the simulations, 
it is found that although P cq (R) is roughly Gaussian and 
(j>(t) is well described by superimposed modes, (R(t)R(0)) 
relaxes much more slowly than the Gaussian expression. 



This is independent of a specific polymer model, relying 
only on an assumption of isotropic Gaussian dynamics. 
The simulated (R(t)R(0)) can be approximately fit for 
long times with Gaussian dynamics by reducing the diffu- 
sion coefficient to a value that also brings the simulations 
into agreement with the bounds on the mean contact time 
provided by the WF and SSS approximations. 

The most fruitful extensions to the Gaussian formal- 
ism would likely include an additional time scale or reac- 
tion coordinate that influences the end-to-end dynamics, 
rather than improving the description of (j)(t) through 
more elaborate Gaussian chain models. Furthermore, 
since the contact rate is directly related to the auto- 
correlation (R(t)R(0)) , it may be more convenient theo- 
retically to focus on this correlation than the intrachain 
mean first passage time to improve polymer dynamics 
models. 



Gaussian Chains 

Although the conclusions (and the analysis of the sim- 
ulations) do not rely on a specific Gaussian polymer 
dynamics model, introducing a rather general Gaussian 
chain at the outset makes the formalism more clear. Con- 
sider a chain consisting of N monomers with positions 
{t*;} and isotropic Gaussian connectivity 

ij 

where (3 = 1/ksT is the inverse temperature, and b is the 
root mean square distance between adjacent monomers 
along the chain. The static correlations of monomer posi- 
tions are determined by the coefficients : (r, -r^/b 2 = 
r^~- , where (• • • ) denotes an equilibrium average over the 
equilibrium distribution Pcq[{^}] ~ cx p[~/3^[{ r }]- (In 
this paper, spatial isotropy of matrices is implicit, e.g., 
Tij is the direct product of an N x N matrix and the 
3x3 indentity matrix.) 

The dynamics are assumed to be over-damped and 
Markovian with a spatially independent friction matrix 
7y. The probability distribution P[{r}, t] evolves accord- 
ing to the Smoluchowski equation 

d t P[{r} ) t}=V[{r}]P[{r},t], (5) 

where the diffusion operator is defined as 

with the diffusion matrix Dij — ksTj^ . 

For a Gaussian chain, the equilibrium correlation func- 
tions and averages for any functions f(R) and g(R) of 
the end-to-end vector R = — f*i can be calculated as 

(/(t)ff(0)) - Jd 3 Rjd 3 R a f(R)g(R )G(Rt\Ro)P cq (Rv) 

(7) 



4 



</) = Jd 3 R.f(R)P cq (R) 7 



(8) 



where the Green's function G(Rt\Ro) and the equilib- 
rium distribution P cq (R) are given by Eq.(Q) and Eq.(||), 
respectively. Note that the chain dynamics are included 
through the dependence of G(Rt\Ro) on the pair corre- 
lation function (j>(t). The form of 4>(t) for this general 
Gaussian chain is given the Appendix. 

End-to-end quenching is incorporated into the Smolu- 
chowski equation as a distance dependent sink term 

d t P[{r},t] = V[{r}]P[{r} 7 t] - ek(R)P[{r} 7 t], (9) 

where ek(R) is the quenching rate for the r.elative end-to- 
end separation R. Following Pastor et aZ.,E3 I consider a 
Heavyside sink 



k(R) 




(10) 



with diffusion limited quenching, e — > oo. The reaction 
rate can be characterized by the mean first passage time 



dtS(t) = S(uj = 0), 



(11) 



where S(ui) is the Laplace transform of the survival prob- 
ability 



S(t) = Jd 3 {r}P[{r},t} 



(12) 



S(lo) = LT t ^[5(t)] = / °°dte- w *S(t). 

The Laplace transformed survival probability (and, 
hence, mean first passage time) for the Gaussian chain 
obey complementary variational bounds. To develop 
this relationship, it is necessary to introduce the adjoint 
Smoluchowski operator, 

ami --£idW£*jw] (13) 

which is defined through T>(P eq f) = —P cq £f, for any 
function f[{r},t}. Writing the probability distribution 
as P(t) = p(t)P cq and assuming equilibrium initial con- 
ditions, Eq.(g) becomes 



[u + £[{r}] + ek{R)]p[{r}M = l 



(14) 



where p — LT t ^ w [p] . In this notation the survival prob- 
ability is S(u>) = (p). 

The variational bounds on the survival probability de- 
rived in Ref. pq can be written in terms of trial functions 
<f[{r}} and £[jr}] 



where the lower bound is given by 

F u [<p] = -{tp(w + £ + ek)ip) + 2(<f), (16) 
and the upper bound is given by@ 
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(?h) + LT t _ e (A m)Ht)] A K(0)fc(0)]) 

(17) 

with A[£(i)fc(i)] = £(t)k(t) - (£fc). The optimal trial 
functions 6F u [<p*] = and 5K U [£*] = satisfy Eq.Q, 
with if* = C = p and F u [<p*] = K u \i*] = S(u). Note 
that Eq. ( |l5| ) evaluated at u> — are bounds on the mean 
first passage time r = 5(0). 

In practice, the bounds on S(ui) and r provided by 
Eqs.([l5|-|l7|) are determined by the choice of trial func- 
tion tOr©ptimize. The upper bound on r was first derived 
by Doic3 who noted that -ftT w= o[£] evaluated at the sim- 
plest trial function £ = Lgives the well known Wilemski- 
Fixman approximation.^!! For the diffusion limited reac- 
tion (e — > oo), the WF approximation becomes 



(18) 



Similarly, the single variable approximation (like SSS) 
can be seen to be a lower bound on r. Restricting the 
trial function to functions of the end-to-end relative vec- 
tor, F u [ip{R)\ can be reduced to averages over P cq (R) 
which are denoted by (• ii)r- This is immediately evi- 
dent for each term in Eq.(|16|) except {ipCip) which reduces 
according to 



/d<p(R) dtp(R) \ 
dri 13 dr.j J 

"J 

/dip(R) d<p{R)\ 
\-^R-- DcS '^R-/ R U)1 

= (<pC(R)<p) n . 



Thus, for this class of trial function the lower bound be- 
comes (for u = 0) 

F u=0 [<p(R)] = -(<p(C{R) + ek)<p) R + 2(ip) R , (20) 
with the single variable adjoint operator 

c w = -i^tt p «w- D "-m> (21) 



and effective diffusion constant 



D eS = D 



NN 



Du-Dut-D 



NX- 



(22) 



Optimizing Eq.(pO|) with respect to ip(R) giveseil 
[C(R) + ek] f = 1. For diffusion limited reactions (e — > 
oo) and a Heavyside sink, the sink term can be rep- 
resented by the boundary condition: <f(R) = 0, for 
\R\ < o,. Thus, the starting_equation introduced by Sz- 
abo, Schulten and Schulteno 



F u [<p]<S{u>)<K u [£], 



(15) 



C(R)f = l, ip(\R\<a)=0, 



(23) 
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is an approximation corresponding to a lower bound for 
the mean first passage time for intrachain contacts. The 
mean first passage time from Eq. (p3|) (r = ^)) can be 
solved exactly giving the SSS approximator 



1 



Tsss 



dx ■ 



1 



D e g J a X 2 P eq (x) 



dy y 2 Pe q (y) 



(24) 



where P eq (R) ~ exp[— /3U(R)] is normalized over a < 
R < oo, and D c s is given by Eq. (B2fl . 

Putting this together, the mean first contact time for 
the original polymer is bounded by 



Tsss < r < t w 



(25) 



Consider applying these bounds to a simulation of a pep- 
tide. Eq. (|25| ) states that the simulated mean contact time 
is bounded by the SSS approximation calculated with the 
simulated potential of mean force and the WF approxi- 
mation calculated from the simulated sink-sink correla- 
tion function. By itself, this is not very helpful in con- 
necting to to microscopic parameters of a polymer model 
that can be generalized to other problems, e.g., poly- 
mer based models of protein folding dynamics. However, 
the bounds calculated from the same many-bead poly- 
mer model are useful since the exact r from the model 
is generally unknown. A chain with Gaussian dynamics 
is chosen in this paper to define the model and associ- 
ated parameters that are compared with the simulated 
contact time and end-to-end dynamics. 

In this case, the two approximations can be simpli- 
fied further. For the SSS approximation, with (3U(R) = 
3R 2 /24-R 2 ), the inner integral of Eq.(p4|) can be evaluated 
givingj23 



T S SS -Doff 



dMu- 3/ VF(3/2,u) 2 , (26) 



6r[3/2,x ] „„, 
where F(a, x) is the incomplete gamma function, and 

3a 2 



X 



2(R 2 



(27) 



In the WF approximation, the sink-sink correlation func- 
tion in Eq. (|l8|) can be calculated using Eq.(|^). For the 
case of a Heavyside sink, Wilemski and Fixman approx- 
imated the pair correlatioxLusing a "unbalanced Heavy- 
side delta approximation'^ that can be written as 



(fc(t)fc(O)) 



(k) 2 erf(y^ ) - (2/ % /7r)yx exp(-a;o)' 

with Xq given by Eq. (p7|) , and 



■CO 



(28) 



(29) 



The approximations Eq.(|26]) and Eqs.(^ 
among those compared with simulations of Gaussian 
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chains by Pastor et a/.;t3 this work verified (though did 
not emphasize) that the two approximations are bounds 
for the mean contact times. 

Since the comparison of these approximations to mea- 
sured (or simulated) contact times is framed in terms 
of the relative value of the diffusion coefficient, this pa- 
rameter warrants closer consideration. As shown in the 
Appendix, a short time approximation to the diffusion 
coefficient also leads to D e ff given by Eq.(|2^). For a free 
draining chain (diagonal friction), D c r = Dn + D\ is 
the relative diffusion constant of the two end monomers, 
which becomes D c s = 2Dq if the friction as uniform 
(D^ = DoSij). Alternatively, Pastor et alr3 also de- 
rived D c g — 2Dq for a Rouse chain (nearest neighbor 
connectivity) with uniform friction using a "local equi- 
librium" approximation, integrating over the internal 
monomer positions in Eq.(|l4|); even though the formal 
manipulations are similar to the variational derivation, 
it is not evident from the "local equilibrium" approxima- 
tion that the SSS approximation with this D e ff is actu- 
ally a lower bound. The diffusion coefficient contained 
the WF approximation follows more naturally from the 
model through the time dependence of the pair correla- 
tion function <f>(t). 

The two approximations to r given by Eq.© and 
Eq.(p4[) suggest (in slightly different ways) that the dy- 
namics of the end-to-end distance R(t) = should 
be give n careful attention. In the SSS approximation, 
Eq. (|23|) describes the dynamics of R(t) in an isotropic 
potential of mean force f3U(R), i.e., R(t) is treated as a 
single reaction coordinate. A single variable treatment 
should be fine to calculate rates as long as the effective 
diffusion constant is appropriate. In the formalism pre- 
sented here, -D ff arises from reducing a more microscopic 
description to an effective one. On the other hand, with- 
out the advantage of a microscopic theory, D c g could 
be determined from the dynamics of the reaction coordi- 
nate (obtained from a simulation, for example). Defined 
in this way, D e g is truly an effective diffusion constant 
of the theory, not necessarily connected to microscopic 
parameters in a straight-forward way (see Ref. 42 for 
example). Turning now to the WF approximation, the 
dynamics of R(t) enters explicitly through the sink-sink 
correlation function, since it is the autocorrelation of a 
function R(t). 

The pair correlation <p(t) has both a distance and an- 
gular component: 4>(t) — (R(t)R(O) cos 9(t))/(R 2 ), where 
cos#(i) is the angle between R(t) and 11(0). Separating 
the distance and angular components, we consider the 
correlation functions 



WW)) 

(R 2 ) ' 



Mt) = (cos0(t)). 



(30) 



The Gaussian Green's function [Eq.(|l])] imposes a specific 
dependence of 4>d{t) and <f>$(t) on <f)(t), just as it does 
on the sink-sink correlation function given in Eq. 
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Performing the integrals in Eq. (]?]), gives 



m = 



l + 2(j> 2 (t) 

m 

2(t) 2 (t) - 1 



arcsin^(i) + 3^1- 4> 2 {t) 



arcsin <f>{t) 



y/l-<P(t) 



(32) 



(The superscript indicates that this is a Gaussian re- 
lationship.) Not only are each correlation a function 
of 4>{t), but these Gaussian expressions obviously also 
have specific relationships to each other, e.g., 3</)|(i) = 

mm)- 

Both <j)\(t) and <pf (t) are normalized to one at t = 0. 
Expanding in powers of <fi(t) gives to lowest order 

1 



(<) 



1 



(33) 



Thus, the Gaussian correlation 4>\{t) approaches the long 
time limit (R) 2 / (R 2 ) = 8/3tt at twice the rate that <j>{t) 
decays, and (pf(t) approaches zero at the same rate as 



Simulation 

In this section, a molecular dynamics simulation of a 
short peptide is presented to investigate the end-to-end 
polymer dynamics and the diffusion limited quenching 
rate approximations developed in the previous section. 
Although simplified models can be very helpful in deter- 
mining the essential underlying physics responsible for 
realistic polymer dynamics, the simulation reported here 
is of a detailed model of the GW 3 peptide studied ex- 
perimentally bjf Lapidus et aln As emphasized by Yeh 
and Hummer J2j simulations of a realistic model provide 
a stringent test for the molecular potentials when com- 
pared directly with measured quenching rates. From a 
more theoretical point of view, simulations obviously also 
make accessible quantities essential to the consistency of 
the theory (e.g., correlation functions) that are not al- 
ways easy to measure in the laboratory. Since this work 
is motivated primarily by the unexplained parameter val- 
ues needed to fit experimentally measured rates, the sim- 
ulation presented in this section is an attempt to capture 
the dynamics relevant to the experimental system. 

The peptide Ace - C(AGQUW - Nme is modeled by 
the force field of Cornell et auB The AMBER 4.1E3 suite 
of programs are used, modified to includeJjhe generalized 
reaction field treatment of electrostaticsoa with a cut- 
off of 9 A. The system, solvated with approximately 1700 
TIP3P water molecules, was held at constant tempera- 
ture (T =r3p0 K) and pressure (P = 1 atm) by Berendsen 
couplings]!!] each with a relaxation time of 0.1 ps. The 
integration step is 0.002 ps. Non-bonded pair list was 
updated every 10 integration steps. A single trajectory 
of 50 ns was generated, starting from an extended con- 
formation, with conformational coordinates saved every 
0.5 ps. 




FIG. 1: Probability distribution of the end-to-end distance 
for CW3. Solid: histogram from simulation with a bin width 
of 0.25 A; dashed: Gaussian distribution with simulated mean 
square distance {R 2 } = 140 A 2 . The vertical dotted lines in- 
dicate the range of contact radii 4 A < a < 6.5 A. 
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FIG. 2: Simulated pair correlation functions of CW3. Solid: 
<f>(t); dotted: cj>e(i); dashed: normalized <j>d (t)r— Rrror bars are 
computed according to Zwanzig and Ai lawadl^l; for example, 
the error bars for <j> e (t) are ±y/2T e /T{\ - 4>$(t)) where T is 
total simulation time (50ns), and Tg = fdt(j)g(t). 



The geometric center of all heavy atoms of the cystine 
and tryptophan sidechains (denoted by cm) are chosen 
to define the end-to-end vector, R(t) = rj m (t) — r° m (t). 
Fig. [l] shows the probability distribution of the end-to- 
end distance R = \R\, as well as a Gaussian distribu- 
tion with the simulated mean square end-to-end distance, 
(■R 2 )sim ~ 140 A 2 . Although there is some deviation, the 
distribution can roughly be described as Gaussian. The 
mean distance (i?) 2 ; m / (R 2 ) S im ~ 0.87 is within a few per- 
cent of the Gaussian relationship (R) 2 / (R 2 ) = 8/3ir 
0.85. 

The simulated pair correlation functions, 4>(t), 4>d(t), 
and <pg(t) are shown in Fig. 0. (Note the radial 
correlation function has been normalized: <t>d(t) = 
(SR(t)SR(0))/(SR 2 ), with SR(t) = R{t)-(R}.) The pair 
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correlation function </>(i) follows the relaxation of the ro- 
tational correlation function <j)g (t) in agreement with the 
limiting Gaussian expression |Eq.(|33|)]. The radial au- 
tocorrelation decays more slowly than the others. Al- 
though a more careful evaluation is considered later, this 
slower relaxation can be seen to be a manifestation of 
non-Gaussian dynamics of the end-to-end relative vector 
R(t) [c.f., Eq.©]. 
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FIG. 3: First passage time distribution of CW3 for a Heavy- 
side sink of radius a = 5 A. The histogram is on a logarithmic 
scale with a bin width equal to powers of 2 (circles); Dashed 
line: exponential distribution, P(t) = AT^ m exp[— r/r s i m ], 
with r sim = {A 2 ) sim /2{A) sim = 2123ps, and A = 0.8. 

The diffusion limited mean first passage time can be 
calculated, from the distribution of time intervals between 
contacts^ {Afc}: 

{A k } = {t s -ti\R(ti-i) <a;R{U)>a;R{tj)<a;}, 

(34) 

for i < j, and k = 1 . . . nj nt . Taking the initial distribu- 
tion to be every conformation with with R(t) > a gives 
the equilibrium distribution of contact times P(t) with 
a mean contact time T s i m = (A 2 ) S j m /2(A) s i m . The simu- 
lated mean first passage times of CW3 for different sink 
radii is given in Table |. For example, a sink radius of 
a = 5.0 A had m n t = 870 separate contact events with a 
mean contact time T s - lm — 2.1 ±0.7 ns. As shown in Fig. [|, 
the distribution of P(r) for a sink radius of a = 5.0 A has 
a strong exponential component (seen by the width of a 
couple of decades) and a non-exponential component at 
short times. A non-exponential component is expected 
because of short-time recrossings at the sharp edge defin- 
ing the Heavyside sink. Although the two components 
are not well separated, the exponential component is 
described reasonably well by P{t) = rT^ cxp(— r/r sim ) 
with an amplitude of 0.8. The simulated distribution 
P(t) for other sink radii are similar. Note that the times 
reported are in the units from the simulation. Anticipat- 
ing the discussion in the next section, to approxmately 
convert this time to units consistent with the viscosity 



plied by a factor of two; for example, the mean time to 
make a contact at a distance a — 4.0 A corresponds to 



7-™* = 5.6 ± 1.4 ns. 



Free Diffusion Coefficient 

A good estimate of the diffusion coefficient of the end 
monomers in TIP3P water is needed to compare the sim- 
ulated contact times with analytical models. To this 
end, two separate 10 ns simulations of H2 — C — OH and 
H2 — W — OH in approximately 1180 TIP3P waters were 
carried out. Otherwise, the conditions are the same as 
those described for the peptide simulation. The value 
determined for the free diffusion coefficient Do appropri- 
ate for the simulation will be used in the next section to 
compare with analytic models. 
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FIG. 4: Simulated diffusion coefficient for cystine (solid) and 
tryptophan (dashed) in TIP3P water. The dotted lines show 
the mean values over 500 ps to 1500 ps, D§ = 0.16 A 2 /ps and 
B% = 0.145 A 2 /ps. 

The translational diffusion coefficient of the geometric 
center of the heavy atoms of the residues can be calcu- 
lated from the Einstein relation D ~ (<5r cm (t) 2 )/6t for 
long times, where 6r cm (t) — r cnl (t) — r cm (0). The diffu- 
sion coefficient for cystine and tryptophan are shown in 
Fig. H Although there is still some variation due to the 
finite time of the simulation, the average diffusion coef- 
ficient over 500ps to 1500ps is taken as an estimate, giv- 
ing D c = 0.16A 2 /ps for cystine and Dq = 0.145A 2 /ps 
for tryptophan. We can verify these values using Stokes 
Law, /cbT / Dq — 67r?7oao, with_ao equal to the van der 
Waals radii of the sidechainsj£3 to find the viscosity of 
TIP3P water, 770. For cystine, using ao — 2.75 A gives 
7/0 = 0.50 cP; and for tryptophan, using ao = 3.38 A gives 
r/o = 0.45 cP; the two estimates differing by 10%. This 
is in reasonable agreement with the viscosity of TIR3P 
calculated by Shen and Freed, ?y = 0.506 ± 0.043 cP.H 

The effective diffusion constant for the single variable 
(SSS) approximation of the mean first contact time is 
taken to be the relative diffusion constant of cystine 
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and tryptophan, D cfi = D c a + = 3.0 x 1(T 5 cm 2 /s. 
Translating this to a uniform monomer friction we have 
D = -Dcff/2 = 1.5 x l(T 5 cm 2 /s. This value for the ef- 
fective diffusion constant is the central assumption in the 
comparison of the simulation with Gaussian chain mod- 
els. 
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TABLE I: Diffusion limited mean first passage times of CW3 
with different sink radii. 



a[A] 


Simulation 

riint a T sim [ns] 6 


-Do = 
tsss [nf 


0.15 A 2 /ps 
] c twf [ns] d 


Do = 

rsss [ns] c 


D /b 
Twf [nsY 


4.0 


816 2.8 ± 0.7 


0.27 


0.67 


1.6 


4.0 


4.5 


901 2.5 ± 0.7 


0.23 


0.61 


1.4 


3.7 


5.0 


870 2.1 ± 0.7 


0.20 


0.56 


1.2 


3.4 


5.5 


1034 1.4 ± 0.3 


0.17 


0.51 


1.0 


3.1 


6.0 


1297 1.3 ± 0.3 


0.15 


0.46 


0.90 


2.8 


6.5 


1495 0.9 ± 0.2 


0.13 


0.42 


0.78 


2.5 



"Number of contact events in simulation. 

^Uncertainty is estimated as the standard error in the mean found 
by the limiting block averaged error using the contact times of each 
configuration, with R(t) > a . 

c Eqs,(|26|-|27|) using the simulated mean square distance, 
(R 2 ) s i^= 140A 2 , and = 2D . 

d Eq. (|18|) and Eqs.(Eq-p3) using the simulated mean square dis- 
tance, (-R 2 ) s i m = 140 A 2 , and pair correlation, </> s im(i)- 

'Same as [c], but with D c g = 2Do- 

■'Same as [d], but with time scaled as t — > t = (Do/£>o)t, i.e., 
using (f> aim (t). 



FIG. 5: Pair correlation <f>{t) for CW3. Lower curves: <f>(t) 
from simulation (solid), and Rouse chain (dotted) with N = 

11, 6 = 3.8 A, and Aj =D Sij,D = 1.5 x 10" 5 cm 2 /s. Upper 
curves, scaled diffusion coefficient, Do = Do/6: <f>(t) from 
simulation with time scaled t — > t = (Do/Do)t (dashed), and 
Rouse chain [Eqs.(p3|-|45[)] with monomer diffusion coefficient 
Do (dotted). Inset is same plot on linear time axis. 

Identifying D e ff as the relative diffusion constant of 
the two end monomers is only rigorously true for free- 
draining chains [Eq.(]22)]. Comparison of the pair corre- 
lation function <j)(t) from simulation and a free-draining 
chain supports for this choice (see Fig. ||). For times 
longer than about 50 ps, <ft B imft) agrees closely with <j>(t) 
from a Rouse chain [Eqs.(|4^-^5|)] with diagonal friction, 
Dij — D Sij. The Rouse chain model uses reason- 
able parameters for this peptide, where the monomers 
of the chain are the a-carbons of the peptide: N = 11, 
b = 3.8 A, and Do = 1.5 x 10~ 5 cm 2 /s is the monomer 
diffusion coefficient determined above. The mean square 
distance of this model Rouse chain, (R 2 ) = (N — l)b 2 = 
144 A, agrees with the simulated value of (-R 2 ) s i m = 
140A. The close agreement of the static and dynamic 
correlations is quite surprising, and one shouldn't ex- 
pect it to hold generally, e.g., for all chain lengths. The 
correspondence with a Rouse chain is considered only 
to support the choice of the effective diffusion constant 
£) e ff = 2D , with D = 1.5 x 10 _5 cm 2 /s — the remain- 
ing analysis assumes only this value for D c ff and Gaussian 
statistics of R(t). 

Except where noted, the reported values for the diffu- 
sion coefficient and contact time correspond to the simu- 
lation viscosity 770 ~ 0.5 cP. To approximately convert to 



the viscosity of water at 298 K and 1 atm, r\ v 



1.0 cP, 



the diffusion coefficient D or contact time r should be 
divided or multiplied by a factor of two, respectively. For 



example, the relative diffusion coefficient corresponding 
to ?7 wat = 1.0 cP is 2D^ at = 1.5 x 10~ 5 cm 2 /s. 



Comparison of Simulation with Gaussian Dynamics 

In this section, the simulated mean contact time and 
correlation functions of CW3 are compared with Gaus- 
sian polymer models. Here, it is demonstrated that the 
effective diffusion coefficient required to approximately 
describe the radial autocorrelation function brings the 
bounding SSS and WF approximations into agreement 
with the simulated mean first contact time. 

The parameters required to calculate the WF and 
SSS approximations can be chosen in a couple of dif- 
ferent ways. For the lower bound, T sss D c g calcu- 
lated with cither the simulat ed p otential of mean force 
U(R) = -kT log Pe q (-R) [Eq.(|lh or the harmonic po- 
tential U(R) = 3R 2 /2(R 2 ) [Eqs.(g§-^)] with (R 2 ) from 
simulation or the Rouse chain agree to within 10% of each 
other for all the sink radii considered. This suggests that 
longitudinal ruggedness in the potential of mean force 
does not significantly reduce the effective diffusion coeffi- 
cient for this peptide^ this conclusion is also consistent 
with the analysis of quenching times forr-thc two shorter 
peptides presented by Yeh and Hummer S3 For the upper 
bound [Eq.(18) with the Gaussian correlations Eqs.(|28|- 
p^)], r WF calculated with either <fi(t) and (R 2 ) from sim- 
ulation or from the Rouse model differ by a few percent. 
Although these calculated rates are in sufficient agree- 
ment for the subsequent analysis, to emphasize indepen- 
dence of a Gaussian chain model, the reported mean first 
contact times use the parameters directly from the simu- 
lation: harmonic potential defined with (i? 2 ) s i m for r sss , 
and 4> s i m (t) and (R 2 ) sim for r WF . 

As shown in Table |, the calculated SSS and WF 
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approximations underestimate the simulated mean first 
contact times for all sink radii considered. The SSS 
approximation is calculated with D e g = 2Dq. If Do 
is treated as a fitting parameter (similar to previous 
analysis of experimentally measured times), the fit to 
r sss defines a new effective monomer diffusion constant 
D = D /a sss with a sss = r sim /r sss . From Table §, 
a sss ranges from approximately 7 for a ~ 6.5 A to 10 
for a — 4.0 A. Following the same reasoning, the WF 
approximation also has an associated effective diffusion 
coefficient. For overdamped motion, an effective diffu- 
sion constant scales time as IDq — tDo, so that us- 
ing Dq as a fitting parameter gives Do = Do/a WF with 
q wf = T s j m /T WP . From Table |, a W F ranges from 2 for 
a = 6.5 A to 4 for a — 4.0 A. (Note that a sss > a WF fol- 
lows from r sss < r WF .) Since SSS and WF are upper and 
lower bound approximations, respectively, the effective 
diffusion coefficient for the original chain should actually 
lie between the values obtained from each approximation; 
then, t sss < T S j m < t wf . For example, for the sink radius 
a = 5.0 A, the D should satisfy D o /10 < D < D /4. 
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FIG. 6: Distance pair correlation (j>d{t) for CW3. Solid: sim- 
ulation; dotted: 4>\{t) using the simulated pair correlation 
</>sim(t) in Eq.(Kll); dashed: same as dotted curve with time 
scaled t — ► i = (Do/Do)t with Do /Do — 6. Inset is same plot 
on linear time axis. 

Why do these approximations require a reduced effec- 
tive diffusion coefficient to fit simulated contact rates? It 
is argued here that the effective diffusion coefficient re- 
flects the slow dynamics of the end-to-end distance. In 
terms of comparing to a Gaussian model, the small effec- 
tive diffusion coefficient accounts for the discrepancy of 
using a Gaussian theory to describe the non-Gaussian dy- 
namics of the end-to-end vector. The simulated correla- 
tions 4>d(t) and 4 > e{^) are shown as solid lines in Fig. ^|and 
Fig. 0, respectively. The dotted lines in these figures show 
the Gaussian expressions 4>\{t) and [Eqs.([3l] p2[ )] 

calculated from the simulated pair correlations 4> S im(t)- 
Comparing the simulated correlations with the Gaussian 
correlations tests the Gaussian statistics of R(t), [i.e., 
this checks whether correlations can be accurately cal- 
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FIG. 7: Angular pair correlation 4>e(t) for CW3. Solid: sim- 
ulation; dotted: (j)^(t) using the simulated pair correlation 
0sim(i) in Eq.(B3); dashed: same as dotted curve with time 
scaled t — > i =(Do/Do)t with Do /Do = 6. Inset is same plot 
on linear time axis. 



culated with the isotropic Gaussian Green's function in 
Eq.(Q)]. Both the simulated radial and angular correla- 
tions decay more slowly than the Gaussian model pre- 
dicts. However, the deviation of the angular correlation 
from Gaussian form is small in comparison to the devia- 
tion of the radial correlations. The semi-log plot in Fig. ^ 
shows that the simulated (/>d(t) and <j)\(t) have roughly 
the same shape, but are translated along the time axis, 
i.e., they differ by a constant scaling of time. As shown 
by the dashed curve in Fig. ||, </>|(t) agrees reasonably 
well with the simulated radial correlation if time scaled 
as t — > t — at with a = 6. For overdamped motion, this 
is equivalent to an effective monomer diffusion coefficient 
Do = Do/ a. As shown in Fig. [| and Fig. |^, this scal- 
ing applied to 4>(t) and <fif(t) results in a poor agreement 
with the simulated correlations. So, rescaling the diffu- 
sion constant to bring the Gaussian and simulated radial 
correlation function into agreement, makes <fi(t) and (j)g(t) 
disagree with the simulated dynamics. 

The dynamics of the end-to-end distance is the most 
relevant to mean first contact times calculations. Ta- 
ble I also shows the calculated times using the effective 
monomer diffusion coefficientEia (Do = Dq/6 with a = 6) 
determined from the fit to 4>d(t). With this scaling, the 
SSS and WF approximations bracket the simulated con- 
tact time, as they should since they are upper and lower 
bounds. Since the approximations provide bounds for a 
rather broad range of a, these results offer only specu- 
lative support that the effective diffusion coefficient de- 
termined by the mean first contact time reflects the dy- 
namics of R(t). Because attention to 4>d(t) arises quite 
naturally from general considerations, the consistency of 
the analysis is very suggestive. 
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Conclusion and Outlook 

The simulated diffusion limited mean first contact time 
agrees with the bounds provided by the WF and SSS ap- 
proximations for a Gaussian chain if the diffusion coef- 
ficient is reduced. When viewed as completely indepen- 
dent approximations, it would seem that the SSS and WF 
approximations would each require different rescalings. 
It was shown, however, that the two approximations are 
complementary variational bounds, giving a range of val- 
ues for the scaling because the approximations should 
merely bracket the measured times, not necessarily fit 
them each independently. 

It was also demonstrated that the relative end-to-end 
distance of the simulated peptide relaxes more slowly 
than the relaxation from Gaussian dynamics. It is ar- 
gued that the effective diffusion coefficient appropriate 
for the calculated mean first contact time also describes 
the end-to-end distance dynamics. Evidence supporting 
this connection is weakened by the broad interval be- 
tween the bounds provided by the SSS and WF approxi- 
mations. Tighter bounds (and other approximations) can 
be generated by the variational formalism through more 
flexible trial functions, helping to make support for this 
relationship more precise. If this relationship is indeed 
valid, one can focus on <pd{t) rather than the mean first 
contact time itself to improve the model. Furthermore, 
the analysis demonstrates that accurately modeling <p(t) 
is not the complete description of the polymer dynam- 
ics, for the Gaussian statistics are checked using <^ s im(t) 
rather the pair correlation from a polymer model. 

It should be noted that Gaussian dynamics and the re- 
duction of a multi-bead model to a single variable model 
are separate issues. The variational bounds r sss < r < 
t wf is a relationship between the mean contact time of 
a single variable model (with a specific D e g) to the cor- 
responding multi-bead model. Reducing the description 
to a single variable model provides the microscopic inter- 
pretation of -D e ff as the relative free diffusion coefficient 
of the end monomers. In general, both models may be 
non-Gaussian, with the potential in the reduced model 
replaced by the potential of mean force of the chain. It 
seems that a microscopic interpretation of the slow end- 
to-end distance relaxation would require a multi-bead de- 
scription, perhaps leading to an effective single variable 
description which is easier to handle analytically. In any 
case, the comparison of the end-to-end distance and an- 
gular relaxation should guide the development of a mi- 
croscopic model that extends beyond the simple Gaussian 
dynamics of a single reaction coordinate. 

In a recent paper, Lapidus et aZ.Lj present further ex- 
perimental work and theoretical modeling of quenching 
in CWfe peptides. In these new experiments, the viscos- 
ity dependence of the quenching rate in CWk peptides 
is measured, isolating the diffusion limited component of 
the measured quenching rate. Fitting the measured diffu- 
sion limited quenching rate with a single variable approx- 
imation corresponding to a stiff chain gives Dq w Dq/10. 



This is the same factor obtained from the fit to the sim- 
ulations presented here with a reaction radius of a < 5 A 
had the SSS approximation been used to define the effec- 
tive monomer diffusion coefficient Dq (see Table |j) . To 
avoid confusion, it is emphasized again that the value 
found from our analysis of the simulation, Dq = -Do/6, 
was determined by the dynamics of the end-to-end dis- 
tance, not the lower bound approximation to the contact 
time. _ 

In the theoretical modeling portion of Lapidus et aZ.,Ej 
the mean first contact time of the end-to-end monomers 
is studied through Langcvin simulations of a model that 
includes excluded volume and chain stiffness (provided 
by the dihedral angle potential), but no other non-local 
interactions between the monomers. The dynamics of 
the squared end-to-end distance from the simulation is 
analyzed in some detail with respect to the single reac- 
tion coordinate description. While the relationship be- 
tween end-to-end distance dynamics and other correla- 
tion functions (angular correlations, for example) is not 
considered, Lapidus et alE3 do emphasize the connection 
between the effective diffusion coefficient appropriate for 
the mean quenching time and the relaxation the end- 
to-end distance. Comparing to experimental data, their 
analysis suggests that the effective diffusion constant de- 
termined by a fit to the single variable approximation is 
not explained by the local properties of the backbone, but 
must involve interactions left out of the Langevin simu- 
lations; e.g., non-local contacts or explicit solvent. This 
very interesting conclusion limits the possible mechanism 
responsible for the slow non-Gaussian relaxation of the 
end-to-end distance. 

The elements missing from these Langevin simulations 
are included in the all-atom molecular dynamics simula- 
tions analyzed in the present paper. Here, an improved 
model with gated diffusion seems promising: the tra- 
jectory of R(t) shows relatively long waiting times sep- 
arating periods of larger diffusional flights in distance 
reminiscent of continuous time random walks, for exam- 
ple. To be consistent with the simulated dynamics, such 
a model should have the property that the gating re- 
tards the end-to-end distance relaxation, but has a rela- 
tively smaller effect on the angular relaxation, at least for 
chains not too long compared to the persistence length. 
More refined analysis of the conformational transitions in 
the simulation should give insight into how to formulate 
a gated diffusion model; this is left for future work. Any 
microscopic model of the gated diffusion should be con- 
strained by the reasonable expectation that the dynamics 
become Gaussian in the limit of very long chains. 
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Appendix: <f>(t) and the short time D c g for Gaussian 
chains 

This appendix gives the pair correlations <fi(t) and de- 
rives the single variable approximation effective diffusion 
coefficient D e ff for a general Gaussian chain. The effec- 
tive diffusion coefficient, derived here as a short time ap- 
proximation, can also be obtained from the variational 
lower bound [Eq.(^2|)] as well as a "local equilibrium" 
approximation^ 

First, the diffusion matrix in the Smoluchowski equa- 
tion [Eq.(||)] is written as a product of a diffusion con- 
stant and a unitless matrix D+j = D^D^ . The Smolu- 
chowski equation can be diagonalized by transforming 
into normal modes rj = ^ p Qi p f p where the transfor- 
mation matrix Q is defined by 



Q _1 6rQ 



Q ^DQ 



-IT 



Q T ro 



diag(Ap) 
diag(fp) 
diag(^p) 



(35) 



with X r 



>p — v p [ip, p = (0 ... -/V — 1). The pair correlations 
of the end to end vector, R = — ri, can then be 
written in terms of the normal modes 



Np ' 



<2i p ) 2 exp 



3£> c 
" b 2 



Apt 



(36) 



Next, consider the effective diffusion constant of the 
single variable approximation to the Gaussian chain. The 
single variable approximation (such as SSS) treats R(t) 
as the only dynamical variable. The Green's function for 
R(t) [Eq.(|])] satisfies the diffusion equation 



d t G(Rt\Ro) = T>(R;t)G{Rt\R ), 
where D is the single variable diffusion operator 







d 



l 



V(R;t)=D(t)-P e<l{ R) dRp ^ {R) 
with a time dependent diffusion coefficient 



D(t) 



3 dt 



log^(i). 



(37) 



(38) 



(39) 



This equation is exact. The single variable approxima- 
tion requires a choice for an effective diffp&ion constant 
D c g to replace D(t). Bicout and SzaboE2 discuss the 
choice of D g and demonstrate some limitations of the 
single variable approximation in the context of rate cal- 
culations of electron transfer in non-Debye solvents. 

The expression for D c g derived from the variational 
lower bound |Eq.(p2|)1 is equivalent to a short time ap- 
proximation, D c fi = D(0). Using A p = v v [i v and 

Aj = J^pVpQipQjp lrom Eq.(|^), 0(0) reduces to 



(R 2 



(40) 



p#0 



3^o 

(R 2 ) 



(D NN + Dn -D 1N -Dm). (41) 



Since 0(0) = 1, D c g = D(Q) becomes 

D e fi = -Daw + fu — Din — D^i 



(42) 



which is Eq.(|2j 

As a specific example, the Rouse chain (nearest neigh- 
bor connectivity) and uniform diagonal friction, Dy — 
DoSij has modes and relaxation rates 



2 - 6, 



N 



P o I p(i ~ 1/2) 
cos 



ttN 



(43) 



X p = V P = 4sm , 
and the mean square end-to-end distance is 

(R 2 ) = (N — l)b 2 . 
The effective short time diffusion constant, 
D cB = 2D , 



(44) 



(45) 



is the relative diffusion constant of the chain ends. Inter- 
nal friction (of the form in Ref. ^||) reduces D c g. 

In the single variable approximation to the end-to-end 
chain dynamics, the pair correlation decays with a single 
time constant 



4>{t) = exp 



3 Jeff 
'{IP) 



(46) 



1 L. J. Lapidus, W. A. Eaton, and J. Hofrichter, Proc. Natl. 
Acad. Sci. USA 97, 7220 (1999). 

2 O. Bieri, J. Wirz, B. Hellrung, M. Schutkowski, 
M. Drewello, and T. Kiefhaber, Proc. Natl. Acad. Sci. USA 
96, 9597 (1999). 



3 R. R. Hudgins, F. Huang, G. Gramlich, and W. M. Nau, 
J. Am. Chem. Soc. 124, 556 (2002). 

4 E. Haas, E. Katchalski-Katzir, and I. Z. Steinberg, 
Biopolymers 17, 11 (1978). 

5 E. Haas, IEEE J. Quantum Electon. 2, 1088 (1996). 



12 



6 M. I. Wallace, L. Ying, S. Balasubramanian, and D. Klen- 
erman, Proc. Natl. Acad. Sci. USA 98, 5584 (2001). 

7 J. D. Bryngelson, J. N. Onuchic, N. D. Socci, and P. G. 
Wolynes, Proteins Struct. Funct. Genet. 21, 167 (1995). 

8 S. E. Jackson, Fold. Des. 3, R81 (1998). 

9 M. Gruebele, Annu. Rev. Phys. Chem. 50, 485 (1999). 

10 L. S. Itzhaki, D. E. Otzen, and A. R. Fersht, J. Mol. Biol. 
254, 260 (1995). 

11 J. J. Portman, S. Takada, and P. G. Wolynes, Phys. Rev. 
Lett. 81, 5237 (1998). 

12 O. V. Galzitskaya and A. V. Finkelstein, Proc. Natl. Acad. 
Sci. USA 96, 112999 (1999). 

13 E. Aim and D. Baker, Proc. Natl. Acad. Sci. USA 96, 
11305 (1999). 

14 V. Munoz and W. A. Eaton, Proc. Natl. Acad. Sci. USA 
96, 11311 (1999). 

15 B. A. Shoemaker, J. Wang, and P. G. Wolynes, J. Mol. 
Biol. 287, 675 (1999). 

16 D. A. Debe and W. A. Goddard III, J. Mol. Biol. 294, 619 
(1999). 

17 J. J. Portman, S. Takada, and P. G. Wolynes, J. Chem. 
Phys. 114, 5082 (2001). 

18 V. Munoz, Curr. Opin. Struct. Biol. 11, 212 (2001). 

19 D. E. Makarov, C. A. Keller, K. W. Plaxco, and H. Metiu, 
Proc. Natl. Acad. Sci. USA 99, 3535 (2002). 

20 R. Zwanzig, J. Chem. Phys. 60, 2717 (1974). 

21 M. Bixon and R. Zwanzig, J. Chem. Phys. 68, 1896 (1978). 

22 C. W. Manke and M. C. Williams, Macromolecules 18, 
2045 (1985). 

23 M. Fixman, J. Chem. Phys. 89, 2442 (1988). 

24 R. Zwanzig, Proc. Natl. Acad. Sci. USA 85, 2029 (1988). 

25 J. D. Bryngelson and P. G. Wolynes, J. Phys. Chem. 93, 
6902 (1989). 

26 D. Thirumalai, V. Ashwin, and J. K. Bhattacharjee, Phys. 
Rev. Lett. 77, 5385 (1996). 

27 S. Takada, J. J. Portman, and P. G. Wolynes, Proc. Natl. 
Acad. Sci. USA 94, 2318 (1997). 

28 T. R. Kirkpatrick and D. Thirumalai, Phys. Rev. B 36, 
5388 (1987). i-i 

29 Very recently) 3 ^ this problem has be reduced to the the 
solution of an integral equation for the special case when 
the quenching has a delta function distance dependence. 

30 I. M. Sokolov, cond-mat/0207159 (2002). 

31 G. Wilemski and M. Fixman, J. Chem. Phys. 60, 878 
(1973). 

32 A. Szabo, K. Schulten, and Z. Schulten, J. Chem. Phys. 
72, 4350 (1980). 

33 R. W. Pastor, R. Zwanzig, and A. Szabo, J. Chem. Phys. 
105, 3878 (1996). 

34 G. Srinivas, A. Yethiraj, and B. Bagchi, J. Chem. Phys. 
114, 9170 (2001). 

35 A. V. Barzykin, K. Seki, and M. Tachiya, J. Chem. Phys. 
117, 1377 (2002). 

36 I.-C. Yeh and G. Hummer, J. Am. Chem. Soc. 124, 6563 
(2002). 

37 L. J. Lapidus, W. A. Eaton, and J. Hofrichter, Phys. Rev. 
Lett. 87, 258101 (2001). 

38 J. J. Portman and P. G. Wolynes, J. Phys. Chem. A 103, 
10602 (1999). 

39 This form of the upper bound is given in Ref. ^8| Eq.(B9) 
with g = £. 

40 M. Doi, Chem. Phys. 11, 107 (1975). 

41 The variational bounds can be similarly developed for the 
single variable equation [C(R) + ek] p — 1. For finite e, 



bounds of this one-dimensional problem generated from 
different trial functions has been addressed in Ref. |38[ We 
do not pursue this approach for Eq.(^3|) because the diffu- 
sion limited rate can be calculated exactly through Eq.(^). 

42 N. D. Socci, J. N. Onuchic, and P. G. Wolynes, J. Chem. 
Phys. 104, 5860 (1996). 

43 W. D. Cornell, P. Cieplak, C. I. Bailey, I. R. Gould, K. M. 
Mertz, Jr., D. M. Ferguson, D. C. Spellmeyer, T. Fox, 
J. W. Caldwell, and P. Kollman, J. Am. Chem. Soc. 117, 
5179 (1995). 

44 D. Pearlman, D. A. Case, J. W. Caldwell, W. S. Ross, 
I. T. E. Cheatam, D. M. Ferguson, U. C. Singh, P. Weiner, 
and P. Kollman, AMBER 4.1 (1995). 

45 L. R. Pratt, G. Hummer, and A. E. Garcia, Biop. Chem. 
51, 147 (1994). 

46 G. Hummer, D. M. Soumpasis, and M. Neuman, J. Phys. 
Cond. Matt. 6, A141 (1994). 

47 H. J. C. Berendsen, J. P. M. Postma, W. F. van Gun- 
steren, A. DiNola, and J. R. Haak, J. Chem. Phys. 81, 
3684 (1984). 

48 R. Zwanzig and N. K. Ailawadi, Phys. Rev. 182, 280 
(1969). 

49 L. J. Lapidus, W. A. Eaton, and J. Hofrichter, J. Mol. Biol. 
319, 19 (2002). 

50 T. E. Creighton, Proteins: Structures and Molecular Prop- 
erties (W. H. Freeman and Company, NY, 1993), 2nd ed. 

51 M. Shen and K. F. Freed, Biophys. J. 82, 1791 (2002). 

52 The agreement, however, is to some degree accidental since 
the SSS approximation depends in part on the weight of 
P cq (R) at small distances, where the Gaussian is a poor 
approximation to the simulated potential of mean force. 
The additional weight of the Gaussian at small R is offset 
to some extent by the peak centered around 4 A in the 
simulated P cq (R). 

53 Converting the effective diffusion coefficient to the viscosity 
77 wat = 1.0 cP, gives 2£>3 rat = 1.5 x 10" 5 cm 2 /s and 2D% at = 
2.5 x 10" 6 cm 2 /s. 

54 L. J. Lapidus, P. J. Steinbach, W. Eaton, A. Szabo, and 
J. Hofrichter, J. Phys. Chem. B (2002), submitted. 

55 D. J. Bicout and A. Szabo, J. Chem. Phys. 109, 2325 
(1998). 



